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Abstract: We investigate non-singular bounce and cyclic cosmological evolutions in a 
universe governed by the extended nonlinear massive gravity, in which the graviton mass 
is promoted to a scalar-field potential. The extra freedom of the theory can lead to certain 
energy conditions violations and drive cyclicity with two different mechanisms: either with 
a suitably chosen scalar-field potential under a given fiducial-metric lapse function, or with 
a suitably chosen fiducial-metric lapse function under a given scalar-field potential. Our 
analysis shows that extended nonlinear massive gravity can alter significantly the evolution 
of the universe at both early and late times. 
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1 Introduction 

The question on whether there exits a consistent covariant theory for massive gravity, 
where the graviton acquires a mass and leads to a modification of General Relativity at 
large distances, was initiated by Fierz and Pauli a long time ago [1]. However, it was 
observed that the nonlinear terms required by massive gravity [2], which can give rise to 
continuity of observables [3, 4], lead inevitably to the existence of the Boulware-Deser (BD) 
ghost [5], and thus make the theory unstable. 

Although for the subsequent decades it was believed that there is no consistent way 
to construct a massive gravity free of ghosts, a nonlinear extension of massive gravity was 
constructed recently by de Rham, Gabadadze and Tolley [6, 7]. In this model, the BD 
ghost can be removed in the decoupling limit to all orders in perturbation theory through 
a systematic construction of a covariant nonlinear action (see [8] for a review). Although 
it is still controversial to prove the absence of BD ghost at the non-perturbative level, 
the theoretical and phenomenological advantages led to a wide investigation of this theory 



Despite the successes of nonlinear massive gravity, it was also noticed that certain cos- 
mological instabilities still exist in the case where the physical and the fiducial metrics have 
simple homogeneous and isotropic forms [42]. This behavior motivated researches towards 
extensions of nonlinear gravity models, namely the construction of nonlinear massive grav- 
ity with less symmetric metrics [13, 43]. However, in [53] a different approach was followed, 
and nonlinear massive gravity was extended allowing for the graviton mass to vary. This 
could be realized by introducing an additional scalar field, which coupling to the graviton 
potentials produces an effective, varying, graviton mass. Moreover, this extension provides 
a natural way to modify General Relativity not only in the IR but also in the UV, since 
the graviton mass can be evolved into a large value at the early universe [54]. Therefore, 
it is interesting to study the cosmological implications of this scenario at early times, and 
this is a main topic of the present work. 



[9-52]. 
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On the other hand, it is well known that cosmological evolution governed by stan- 
dard Einstein gravity usually suffers from the problem of initial singularity if Null Energy 
Condition (NEC) is preserved [55]. A potential solution to the cosmological singularity 
problem may be provided by non-singular bouncing cosmologies [56-58]. Such scenarios 
have been constructed within various approaches to modified gravity [59-73]. Alternatively, 
a non-singular bouncing model can be acquired by using NEC violating matter [74-78], 
by making use of various mechanisms [79-86]. Furthermore, the perturbation theory of 
non-singular bounce cosmology and its relation to observables was extensively studied in 
the literature [76-78, 87-92]. The extension of all the above bouncing scenarios is the old 
idea of cyclic cosmology [93], in which the universe presents a periodic sequence of contrac- 
tions and expansions. Cyclic cosmology has attracted a significant interest the last years 
[94] since it brings different insights for the origin of the observable universe [95-106] (see 
[107, 108] for recent reviews). 

In the present work we are interested in constructing scenarios of cyclic cosmology in a 
universe governed by the extended, varying-mass, nonlinear massive gravity. Particularly, 
we first determine the fiducial-metric lapse function and we suitably reconstruct the form 
of the potential of the scalar field which leads to a cyclic universe. Alternatively, we first 
determine the scalar potential and we reconstruct the fiducial-metric lapse function in order 
to obtain cyclicity. Interestingly enough cyclicity is easily acquired in this framework, since 
extended nonlinear massive gravity can violate certain energy conditions and therefore has 
fruitful implications to physics of the universe at both early and late times. 

This paper is organized as follows. In Section 2 we briefly review the cosmological 
equations under the extended scenario of nonlinear massive gravity with a scalar field being 
introduced to evolve the graviton mass. In section 3 we construct scenarios of bouncing 
and cyclic universe. In particular, in subsection 3.1 we start with a given fiducial-metric 
lapse function and we reconstruct the scalar potential; while in 3.2 we start from a given 
scalar potential and we determine the corresponding fiducial-metric lapse function that 
leads to cyclicity. Finally, section 4 is devoted to the summary of our results. 



2 Cosmology in extended nonlinear massive gravity 

In this section we briefly review cosmology in extended nonlinear massive gravity [53, 54]. 
In such a scenario the graviton mass is upgraded into an evolving function depending on a 
cosmic scalar field, of which a canonical kinetic term and a standard potential are added 
in the action too. In this construction, the total action is written as 



S 



+ V(^)(U 2 + a 3 U 3 + a 4 U A ) - \d^d^ - W[$) 



+ S m , (2.1) 



where M p is the reduced Planck mass, R is the Ricci scalar, and ip is the newly introduced 
scalar field with W{ip) its usual potential and V(ip) an additional potential coupling to the 
graviton potentials. These graviton potentials are given by, 



U2 = ic^K) , 



u 4 



(2.2) 
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with 

n = K~ yJg^fABdp^d^ , IC^^^K^-IC^), (2.3) 

and similarly for the other antisymmetric expressions. Moreover, 03 and 04 are dimen- 
sionless parameters. Additionally, /ab is a fiducial metric, and the four-form fields 4> A (x) 
are the Stiickelberg scalars introduced to restore general covariance [109]. As it was shown 
in [53] the above scenario is still free of the the BD ghost. Finally, in the action (2.1) we 
can take into account the standard matter action S m , minimally-coupled to the dynamical 
metric, corresponding to energy density p m and pressure p m . 

In order to derive explicitly the cosmological equations, we need to impose certain 
ansatzes for the two metrics. For the physical metric we consider a flat Friedmann- 
Robertson- Walker (FRW) form (one can also straightforwardly investigate the non-flat 
case [53]): 

d 2 s = -N(t) 2 (It 2 + afrfSijdafda?, (2.4) 

with N(t) the lapse function and a(r) the scale factor, and for simplicity for the Stiickelberg 
fields we choose the forms as follows, 

0° = b(r), ft = x\ (2.5) 

Note that although in standard massive gravity such a choice for the dynamical metric 
cannot be accompanied by a simple choice for the fiducial one, in the scenario at hand the 
extra freedom allows for a simple Minkowski ansatz for the fiducial metric 

Iab = VAB- (2.6) 

Variations of the action with respect to N and a give rise to the following Friedmann 
equations 

3M 2 P H 2 = PMG + Pm, (2.7) 
-2M 2 P H = p MG + PMG + Pm+ Pm, (2.8) 

where we have defined the Hubble parameter H = a/ a, with a = da/(Ndr). In the end we 
adopt N = 1 for convenience. In the above expressions we have defined the energy density 
and pressure arising from the modified gravitational sector as 1 

PMG = \i? + WW + V{i>) Q - l) lh(a) + h(a)} (2.9) 

PMG = \i? ~ Wty) ~ V(i>)h(a) - W)6/i(a) , (2.10) 

1 Note that there exists a \ factor in the expressions of energy density and pressure as well as in the 
functions fi, /i, fs, and f±. The value 1 is the default reference that we have chosen in the fixed physical 
metric of form (2.4) and the Stiickelberg fields of form (2.5). One can rescale the reference value as an 
arbitrary one a re /, by replacing a by a/a re f in the metric or transforming (j? to 4> l /a re f in the Stiickelberg 
fields, but for simplicity throughout this work we set a re f = 1. 
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In addition, the usual continuity equation is still preserved: 

Pmg + 3H(p M G + Pmg) = 0. 



(2.11) 



(2.12) 



Variation of the action with respect to the scalar field ip provides its evolution equation: 



* + ^ + ^ + ^{(H-i)lAW + A(«)] + 3iAW 

Additionally, variation with respect to b provides the constraint equation 

V(iP)Hf 1 (a)+V(iP)f 2 (a) = 0, 
which, as it was shown in [54], using (2.11) in general leads to 

V(^(T)) = V e 



0. 



f -^da Vq 

J a J2 = 

(a — 1)[«4 — (3q3 + 2a4)a + (3 + 3Q3 + a^a 2 ] 



(2.13) 



(2.14) 



(2.15) 



This relation restricts radically the coupling-potential V(tp) 2 and can be very helpful since 
it offers the behavior of V(^(r)) without the need to find explicitly the solution ij}{r) [which 
obviously is consistent with this V(i/j(t))]. One could still worry about the fact that V(ip) 
becomes negative for a < 1 (or a < a re f in general, which can be push to very small values 
but still non-zero) , however the bouncing behavior of the present work offers a solution to 
this problem, since setting the bounce scale factor to be larger than 1 (or a re f in general) 
ensures that the graviton mass square will be always positive. Finally, the equations close 
by considering the matter evolution equation p m + 3H(p m + p m ) = 0. 



3 Bouncing and Cyclic solutions 

The above cosmological scenario proves to exhibit very interesting behavior. In particular, 
the form of the coupling-potential (2.15) implies that in an expanding universe V(iP(t)) 
always goes to zero at late times, and thus the scenario at hand always gives the stan- 
dard quintessence scenario [54]. The only case that this will not happen is if the scalar 
field dynamics is so effective that it will change the universe evolution from expansion to 

2 Note that in standard massive gravity this equation imposes strong constraints on b(r) [13, 16, 51], but 
in the present construction the extra freedom brings these constraints on V(ip), leaving b(r) free [54]. 
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contraction. In this case V{iP{t)) will start increasing and at some point it can trigger a 
nonsingular bounce. A successive sequence of bounces and turnarounds offers the cyclic 
cosmology. In this section we analytically explore these possibilities. 

Whether a universe is expanding or contracting depends on the positivity of the Hubble 
parameter, that is, in the contracting phase H is negative while in the expanding one it 
is positive. By making use of the continuity equations it follows that at the bounce and 
turnaround points H = 0, however at and around the bounce we have H > 0, while at and 
around the turnaround we have H < 0. If a bounce and a turnaround can be realized in a 
given cosmological scenario, then the aforementioned general requirements for H must be 



Interestingly enough, observing the form of the two Priedmann equations (2.7), (2.8) 
along with (2.9), (2.10), (2.15), we can see that the above requirements can be easily 
fulfilled in the present scenario, for suitable potentials W{ip) and/or suitable forms for b(r). 
In particular, the coupling-potential V{ip) can trigger the bounce and the turnaround, and 
actually it is this term that can cause the necessary null energy condition violation, as it 
was mentioned in [54], where this term could lead the effective dark energy caused by the 
gravitational modification to exhibit phantom behavior. 

In summary, cosmology governed by extended massive gravity has a large freedom 
to fulfill the above bounce and cyclicity requirements. In particular, for a given 6(r) 
we can adjust the potential W{ip) in order to obtain a desired scale-factor evolution, or 
alternatively for a given W{ip) we can adjust b(r). In the following subsections we examine 
these two cases separately. 

3.1 Known b(r) 

Suppose we determine 6(r) at will. Let us first start from the desired result, that is we 
impose a known form of the scale factor a(r) possessing a bouncing or oscillatory behavior. 
In this case both H{t) and H(t) are straightforwardly known. Therefore, we can use the 
Friedmann equations (2.7), (2.8) together with (2. 9), (2. 10), in order to extract the relations 
for tp(r) (through ?P(t)) and W{t), obtaining 



where p m (a(r)) and p(a(r)) are known as long as the matter equation-of-state parameter 
w m = Pm/p m is determined. Then, eliminating time between these two expressions we can 
extract the explicit form of the potential W(ifj), which is the one that generates the desired 
a(r)-form. Finally, we mention that the basic requirement for the above procedure to be 
valid and a bounce or cyclic behavior to be possible, is that the a(r) and 6(r) ansatzes and 
the parameter choices must lead to a positive V> 2 (t). 



fulfilled. 




(3.1) 



W{t) = M 2 p 3H\t) + 2H(t) + Pm {a{T)) 



V(^(t)) h{a{r)) + b{T)h{a{T)) 



(3.2) 
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We now proceed to a specific, simple, but quite general example. Firstly, we consider 
6(r) to have the simple form 

6(r) = r, (3.3) 

and thus b = 1 in the above expressions. Now, in order to get an idea of what kinds of 
potentials W{tp) could lead to cyclic behavior, we perform a test-procedure starting from 
the desired result, that is we assume a cyclic universe with an oscillatory scale factor of 
the form 

a(r) = Asm(ojT) + a c , (3-4) 

where the non-zero constant a c is inserted in order to avoid the singularity. Thus, r 
varies between — oo and +oo, with r = being a specific moment without any particular 
physical meaning. Finally, note that the bounce occurs at as = a c — A. From (3.4) we 
straightforwardly find H{t) and H(t), and thus substitution into (2.15) and then to (3.1) 
and (3.2) gives the corresponding expressions for ip(r) and W(r). 

In order to provide a more transparent picture of the obtained cosmological behavior, in 
Fig. 1 we present the evolution of i/j(t) and W(t) for the scale factor (3.4) with A = 5, u = 1 
and a c = 7, assuming also dust matter (w m = 0), that is p m = and p m {a) = /O m B«s/ a3 ' 
with p m B the energy density at the bounce. Thus, eliminating time between tp(r) and W{t) 

10 

> 
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Figure 1. The evolution o/i/'(r) and W(t) as functions of time in the case where b(r) = r and 
a(r) = As'm(ujT) + a c . with A = 5, uj = 1, and a c = 7. In numerical elaboration we use = 1. 
CK4 = 1, Vq = 11, PmB = 0.01, and Mp = 0.1. All dimensional parameters are normalized in unit 
ofM p . 

allows us to re-construct the corresponding expression for W(ijj), shown in Fig. 2. From 
these figures we observe that an oscillating and non-singular scale factor can be generated 
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Figure 2. The reconstructed potential W(ip) in the case where b(r) = r and a(r) = Asui(lot) + a c , 
with A = 5, lo = 1, and a c = 7. In numerical elaboration we use 013 = 1, on = 1, Vo = 11, 
PmB = 0.01, and Alp = 0.1. All dimensional parameters are normalized in unit of M p . 



by an oscillatory scalar potential W(ip). This W(^>)-form was expected, since a non- 
oscillatory W(ijj) would be difficult to produce an infinitely oscillating scale factor. Finally, 
we mention that although we have presented a specific example, we can straightforwardly 
repeat the described procedure imposing an arbitrary oscillating form for the scale factor. 

The aforementioned bottom-to-top approach was enlightening about the form of the 
scalar potential that leads to a cyclic cosmological behavior. Therefore, one can perform 
the above procedure restoring cause and effect, that is starting from a specific oscillatory 
W(ip) and resulting to an oscillatory a(r). In particular, (3.1) is written in a compact 
form as ip 2 {r) = Pi(a,a,a) and similarly (3.2) as W(t) = P2(a,a,a). Thus, we can 
invert the known form of W(ip) = W(iP(t)) acquiring i/j(t) = W^~^ (P2{a, a, a)), and then 
iP 2 ( T ) = {£ [W^ (P 2 (a,a,a))]} . Therefore, the scale factor arises as a solution of the 
differential equation 



Pi (a, d, d) 



d 



W±-V (P 2 (a,a,a)) 



(3.5) 



As a specific example we consider the simple form 



W(iJ>) = W sm(co w -0) + W c 



(3.6) 



In this case ijj 



sm 



W(i/j(t))-W c 
Wo 



where W(iJj(t)) = W{t) = P 2 (a,a,a) with 



P2(a,a,a) the right hand side of expression (3.2). Therefore, the differentiation leads 
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to: 



V(r) 



1 



Wqujw 



— [P 2 (a,d,d)] 



P 2 {a,a,d)-W c 
W 



2 dr 



(3.7) 



and thus we obtain 



Pi (a, d, a) = < 



1 



Wq LOW 



1 



P 2 (a,d,a)-W e 
Wo 



2 dr 



[P 2 (a,d,d 




(3.8) 



where as we have mentioned Pi (a, d, d) is the expression inside the curly brackets in (3.1). 

Differential equation (3.8) cannot be solved analytically, but it can be easily elaborated 
numerically, and moreover inserting its solution for tt(r) into (3.7) we obtain VK 7 ") too. In 
Fig. 3 we depict the corresponding solutions for a(r) and ip( T )i f° r the ansatz (3.6) with 
Wo = 0.2, u w = 3 and W c = 0.2, for a 3 = 1, a 4 = 1, Vb = 11, M P = 0.1, p mB = 0.01, in 
M|,-units (the potential parameters have been chosen in order to acquire a cyclic universe 
with as ~ 2 at the bounce, similarly to the previous example). In conclusion, we verified 
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Figure 3. TTie evolution of the scalar field t/j(t) (upper graph) and of the scale factor (lower graph), 
for b(r) = t and for the scalar-field potential W(ij)) = Wq sin(wvv ip) + W c with Wq = 0.2, ujw = 3 
and W c = 0.2. In numerical elaboration we use = 1, = 1, Vq = 11, p m B = 0.01, and 
Alp = 0.1. All dimensional parameters are normalized in unit of M p . 



that an oscillatory scalar-field potential can lead to a cyclic universe. 

We close this subsection by mentioning that here b(r) was considered by hand and 
then we determined the potential W{ip) in order to obtain cyclic behavior. The only 
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requirement for the above analysis to be valid is to obtain a non-negative ip 2 in (3.1). If 
this is not possible we deduce that the corresponding parameter-choice and /or the imposed 
6(r)-ansatz cannot lead to a cyclic solution. However, note that there are 6(r) ansatzes 
that cannot lead to cyclicity independently of the parameter choice (for a3,a4,Vo,p m B,w m ). 
In particular, using the forms 6(r) = a(r) and 6(r) = In a(r) (which in [54] it was shown 
to lead to phantom behavior) we could not find any parameter choice that could lead to 
non-negative ip 2 and therefore to a cyclic behavior. Thus, cyclicity cannot be obtained for 
every 6(r)-form. 

3.2 Known W(ip) 

Suppose we determine W(ip) at will, and similarly to the previous subsection as a first 
test-procedure we impose a known scale factor o(r) possessing oscillatory behavior, that 
is both H(t) and H(r) are known. In this case eliminating ip between the two Friedmann 
equations (2. 7), (2. 8) straightforwardly gives 

1 



b(r) = — +[V(a(T))f 1 (r)}- 1 \2M P H(t) + 3H z (t) - p m (a(r)) + p m (a(r)) 



o(r) 

-2W{iI>{t)) - 2V(a(r)) 



2/ 



1 

1 



a(r) 



[/3(a(r)) + A(a(r))] = P 3 (^r), (3.9) 



where V(cl{t)) is given by (2.15), and thus substitution to the scalar-field evolution equation 
(2.13) 
form 



(2.13), along with the chain rule ^ = leads to a differential equation for ip of the 



i> + 3Hj> + ^ + ~ | Q - l) [/ 3 (a) + /i(a)] + 3P 3 (t/>, r)/ 2 (a)| ss P 4 (V, ^, <) = 0, 

(3.10) 

that can be numerically solved, leading to V'( T )- Thus, inserting ^(r) back to (3.9) gives a 
simple differential equation of the form b = P§{t) that can be easily solved as 



b(r) = J P 5 (r')dr'. (3.11) 

As a specific example we consider the simplest case in which W(ip) is constant, namely 

W(ip) = W . (3.12) 

Similarly to the previous subsection we assume a cyclic universe with the oscillatory scale 
factor (3.4), namely a(r) = Asui(ujt) + a c . In this case from (3.9) we can see that b(r) is 
analytically known, which simplifies the procedure. Assuming also dust matter (w m = 0) 
and p m (a) = p m B<^Bl a ^^ with p m B the energy density at the bounce, we numerically 
solve the differential equation (3.10), in order to extract tp(r). In Fig. 5 we depict the 
corresponding evolution of ip(r), b(r) and b(r), and thus we deduce that this b{r) form 
generates the oscillatory universe. As we observe, 6(r) has an oscillatory profile, although 
not of a simple form, being always positive, and this leads b(r) to be monotonic, but with 
varying growth rate. The fact that b(r) and not b{r) should be oscillatory in order to 
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Figure 4. The evolution of the scalar field ip(r) (upper graph) , of the time- derivative of the fiducial- 
metric lapse function b(r) (middle graph) and of the fiducial-metric lapse function 6(r) (lower 
graph), in the case where W(ip) = Wq with Wq = 0.02. In numerical elaboration we use = 10, 
«4 = 10, Vo = 500, p m B = 0.05, and Mp = 0.1. All dimensional parameters are normalized in unit 
ofM p . 



produce an oscillatory scale factor was expected, since in the cosmological equations it is 
b and not b that appears. 

We stress here that the differential equation (3.10) in the above procedure cannot be 
solved for every W(ip), and additionally one has to suitably choose the model parameters. 
Thus, although in the previous subsection one had full control on which cases lead to 
oscillatory solutions (namely those that give a positive ip 2 ( T ) m (3-1)), in the present 
subsection one cannot easily see which W(tp) ansatzes and/or which parameter choices can 
lead to cyclic solutions. 

Finally, since from the above bottom-to-top approach we deduced that in order to get 
an oscillatory universe we must use an oscillatory 6(r), in the following we restore cause 
and effect and we impose such a 6(r) as an input, extracting the resulting ip(r) and a(r) 
from the solutions of the cosmological equations. Although this procedure is possible in 
general, in the following we apply it in the simple example of a constant W(ip) = Wo of 
(3.12). In this case, and with a given b(r), equation (3.9) is a simple differential equation 
for a(r), which can be numerically solved. Then, inserting this a(r) back in (3.1), we find 
also the solution for ip(r). Similarly to the previous subsection, the only requirement for 
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this procedure to work is to obtain a positive ip 2 (r) . 

As a specific example for b(r) we impose the simple ansatz 

6(r) = B sm(u B T) + b c , (3.13) 

that is we impose a fiducial-metric lapse function of the form 

b(r) = cos (cj b t) + b c r + c, (3.14) 

lob 

with c an integration constant. We then solve numerically (3.9) for a(r) and (3.1) for ?/>(r) 

and in Fig. 5 we depict the corresponding results. Indeed, the fiducial-metric lapse function 




t 1 1 • 1 1 1 • r 



-1.0 -0.5 0.0 0.5 1.0 

4n 




2 i 1 1 1 1 1 1 1 r~ 

-1.0 -0.5 0.0 0.5 1.0 

X 

Figure 5. The evolution of the scalar field ip(r) (upper graph) and of the scale factor (lower graph), 
in the case where W(ip) = Wo and b{r) = B sm(uJBT) + b c , with Wo = 0.001, B = 0.1, u>b = 0.01 and 
b c = 2.3. In numerical elaboration we use — 4, a± — 4, Vo = 500, p m B = 0.05, and Mp = 0.1. 
All dimensional parameters are normalized in unit of M p . 

(3.14) generates an oscillatory universe, although not of a simple form since the turnarounds 
are softer than the bounces. However, we mention that one needs to carefully select the 
model parameters in order to obtain an oscillatory scale factor under the 6(r)-ansatz (3.14), 
however he has control of the procedure, since he needs to fulfill the requirement of V' 2 (r) 
to be positive. 
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4 Discussion 



In the present work we investigated bouncing and cyclic cosmological behaviors in a uni- 
verse governed by extended massive gravity, in which the graviton mass has been promoted 
to a function of an extra scalar field. This model has additional freedom comparing to usual 
massive gravity, and thus it leads to significantly different and richer behavior [53, 54]. In 
particular, although the scenario at late times tends to coincide with standard quintessence, 
at early and intermediate times the effective graviton mass can be large and thus play a 
crucial role in the universe evolution. Amongst others, the capability of the scenario to 
lead to violation of the Null Energy Condition can lead to a bounce or a turnaround, the 
successive sequence of which can naturally give rise to cyclic cosmology. 

Extended nonlinear massive gravity has an enhanced freedom (apart from the extra 
scalar field one can see that the fiducial-metric lapse function b(r) is not constrained as in 
usual massive gravity and it can be free), therefore it can drive cyclicity with two different 
mechanisms. The first is to impose an arbitrary fiducial-metric lapse function 6(r) and 
suitably choose the usual scalar field potential W(ip) in order to obtain a cyclic behavior. 
It proves that one should use an oscillating W{ip) as expected, and the only requirement 
for this procedure to hold is to obtain a positive ?/> 2 (t). Thus, this is not possible for every 
b(r) form, that is not every b(r) can be consistent with cyclicity. 

The second mechanism to drive cyclic behavior is exactly the fiducial-metric lapse 
function b(r). In particular, imposing an arbitrary scalar-field potential W(ip) we suitably 
choose 6(r) in order to obtain a cyclic behavior. It proves that one should use an oscillating 
6(r) as expected (since 6(r) appears in the equations and not b(r)). Similarly to the first 
mechanism above, not all scalar potentials are consistent with cyclicity. 

A crucial issue in all bouncing and cyclic scenarios is the processing of perturbations 
through the bounces. A simple check on the stability of this type of models is to look at 
the generalized Higuchi bound derived in [48] (see also [110] for a different expression in an 
earlier model). Such an analysis in the scenario at hand shows that the effective graviton 
mass nearby the bounce is characterized by the parameter Vo, which is approximately 
between O(10) and 0(100) in M 2 units, that is much higher than the generalized Higuchi 
bound. Therefore, at this level, we can conclude that our model is well behaved when the 
perturbations pass through the bouncing points. However, we should notice that since a 
cosmic scalar is introduced to drive the graviton mass varying along background evolution, 
the stability issue arisen from this scalar field ought to be taken into account in a global 
analysis. Such a complete perturbation analysis of the extended nonlinear massive gravity 
lies beyond the scope of the present work and it is left for future investigation. 
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